from scipy import ndimage
from skimage import data,util
from matplotlib import  pyplot as plt

img=data.camera()
#统计排序滤波
#对图像加入椒盐噪声
noise_img=util.random_noise(img,mode='s&p',seed=None,clip=True)
#中值滤波
n=3
new_img=ndimage.median_filter(noise_img,(n,n))
#最大值滤波
max_img=ndimage.maximum_filter(noise_img,(n,n))
#最小值滤波
min_img=ndimage.minimum_filter(noise_img,(n,n))

#显示图像
plt.figure()
plt.imshow(img,cmap='gray')#原始图像
plt.figure()
plt.imshow(noise_img,cmap='gray')#加噪处理图像
plt.figure()
plt.imshow(new_img,cmap='gray')#中值滤波处理图像
plt.figure()
plt.imshow(max_img,cmap='gray')#最大值滤波处理图像
plt.figure()
plt.imshow(min_img,cmap='gray')#最小值滤波处理图像



#空间滤波--平滑处理
#定义二维灰度图像的空间滤波函数
def correl2d(img,window):
    #使用滤波器实现图像的空间相关
    #mode='same'表示输出尺寸等于输入尺寸
    #boundary='fill'表示滤波前，用常量值填充原始图像的边缘，默认常量值为0
    s=signal.correlate2d(img,window,mode='same',boundary='fill')
    window1=np.ones((3,3))/(3**2)#3*3盒状滤波模板
    window2=np.ones((5,5))/(5**2)#5*5盒状滤波模板
    window3=np.ones((9,9))/(9**2)#9*9盒状滤波模板
    return s.astype(np.uint8)

new_img1=correl2d(img,window1)
new_img2=correl2d(img,window2)
new_img3=correl2d(img,window3)
    
plt.figure()
plt.imshow(img,cmap='gray')#原始图像
plt.figure()
plt.imshow(new_img1,cmap='gray')#3*3盒状滤波图像
plt.figure()
plt.imshow(new_img2,cmap='gray')#5*5盒状滤波图像
plt.figure()
plt.imshow(new_img3,cmap='gray')#9*9盒状滤波图像



#高斯平滑滤波
def correl2d(img,window):
    #使用滤波器实现图像的空间相关
    #mode='same'表示输出尺寸等于输入尺寸
    #boundary='fill'表示滤波前，用常量值填充原始图像的边缘，默认常量值为0
    s=signal.correlate2d(img,window,mode='same',boundary='fill')
    return s.astype(np.uint8)
#    img=data.camera()
    #定义二维高斯函数
    def gauss(i,j,sigma):
        return 1/(2*math.pi*sigma**2)*math.exp(-(i**2+j**2)/(2*sigma**2))
    #定义radius*radius的高斯平滑函数
    def gauss_window (radius,sigma):
        window=np.zeros((radius*2+1,radius*2+1))
        for i in range(-radius,radius+1):
            for j in range(-radius,radius+1):
                window[i+radius][j+radius]=gauss(i,j,sigma)
        return window/np.sum(window)
    window1=gauss_window(3,1.0)#3*3高斯平滑滤波模板
    window2=gauss_window(5,1.0)#5*5高斯平滑滤波模板
    window3=gauss_window(9,1.0)#9*9高斯平滑滤波模板
#生成滤波结果
new_img1=correl2d(img,window1)
new_img2=correl2d(img,window2)
new_img3=correl2d(img,window3)
#显示图像
plt.figure()
plt.imshow(img,cmap='gray')#原始图像
plt.figure()
plt.imshow(new_img1,cmap='gray')#3*3高斯平滑滤波图像
plt.figure()
plt.imshow(new_img2,cmap='gray')#5*5高斯平滑滤波图像
plt.figure()
plt.imshow(new_img3,cmap='gray')#9*9高斯平滑滤波图像


#拉普拉斯锐化

def laplace_sharpen(image,c):
    '''
    拉普拉斯锐化
    ：param image:输入图像
    ：return :输出图像
    ：param c:锐化系数
    '''
    
    input_image_cp=np.copy(image)#输入图像的副本
   #拉普拉斯滤波器
    laplace_fliter=np.array([
            [1,1,1],
            [1,-8,1],
            [1,1,1],
            ])
    input_image_cp=np.pad(input_image_cp,(1,1),mode='constant')#填入输出图像
    m,n=input_image_cp.shape#填充后的图像的尺寸
    output_image=np.copy(input_image_cp)#输出图像
    for i in range(1,m-1):
        for j in range(1,n-1):
            R=np.sum(laplace_filter*input_image_cp[i-1:i+2,j-1:j+2])#拉普拉斯滤波器响应
            output_image[i,j]=input_image_cp[i,j]+c*R
    output_image=output_image[1:m-1,1:n-1]#裁剪
    
    return out_image

#图像空间滤波函数
def correl2d(img,window):
    
    m=window.shape[0]
    n=window.shape[1]
    #边界通过0灰度值填充扩展
    img1=np.zeros((img.shape[0]+m-1,img.shape[1]+n-1))
    img1[(m-1)//2:(img.shape[0]+(m-1)//2),(n-1)//2:(img.shape[1]+(n-1)//2)]=img
    img2=np.zeros(img.shape)
    for i in range(img2.shape[0]):
         for j in range(img2.shape[1]):
             temp=img1[i:i+m,j:j+n]
             img2[i,j]=np.sum(np.multiply(temp,window))
    return img2

#img为原始图像
imge=data.camera()
#img_laplace为原始图像经过拉普拉斯变换后的结果
window=np.array([[-1,-1,-1],[-1,8,-1],[-1,-1,-1]])
img_laplace=correl2d(imge,window)
img_laplace=255*(img_laplace-img_laplace.min())/(img_laplace.max()-img_laplace.min())
#将img和img_laplace相加得到锐化增强图像img_laplace
img_laplace_enhance=imge+img_laplace
#img_sobel为原始图像img进行sobel处理的结果
img_sobel=filters.sobel(imge)
#使用5*5均值滤波器平滑后的sobel图像
window_mean=np.ones((5,5))/(5**2)
img_sobel_mean=correl2d(img_sobel,window_mean)
#将img_laplace_enhance*img_sobel_mean相乘得到锐化结果
img_mask=img_laplace_enhance*img_sobel_mean
#将原始图像img 和锐化结果img_sharp相加得到锐化增强图像
img_sharp_enhance=imge+img_mask
#对img_sharp_enhance进行灰度幂律变换得到最终结果
img_enhance=img_sharp_enhance**0.5
#显示图像
imgList=[imge,img_laplace,img_laplace_enhance,img_sobel,img_sobel_mean,img_mask,img_sharp_enhance,img_enhance]
for grayImg in imgList:
    plt.figure()
    plt.axis('off')
    plt.imshow(grayImg,cmap='gray')
plt.imshow(img_enhance,cmap='gray')